Multiple Alignment with Hidden Markov Models

Multiple alignment of sequences is an important problem in bioinformatics. For example, multiple alignment of proteins belonging to the same family can provide valuable information for the common protein structure. Looking at a multiple alignment one can identify patterns which are not obvious if one looks at just pairwise alignments between members of the family. Traditional alignment schemes first choose a scoring matrix and gap penalties and then use dynamic programming which gives exact result in O(N) operations for two sequences of length N. For multiple alignment of K sequences dynamic programming requires O(N) operations, which is exponential in K. Therefore the algorithm is not feasible for multiple alignment of more than a few sequences. On the other hand, protein families often contain hundreds of sequences. To avoid the problem various heuristic algorithms were introduced. However, it is often the case that the solution to a multiple alignment problem depends strongly on the scoring scheme. This poses the question of how to choose the most biologically relevant one. Another very important issue is that the scoring matrices and gap penalties are assumed to be independent of position along the sequence, while gene families often exhibit regions which are highly conserved and regions which are highly variable. In order to obtain the most biologically relevant multiple alignment these differences should be taken into account. If we want to construct multiple alignments based only on primary structure information, the variations of scoring matrices and gap penalties should be deduced from the sequences themselves. Hidden Markov Models (HMMs) are an implementation of the idea that the scoring parameters should guide the multiple alignment as much as the alignment should determine the scoring parameters. A Hidden Markov Model (HMM), λ, is a stochastic machine for generating (amino acid) sequences. Different sequences are generated with different probability. A sequence S is generated with probability P (S|λ). The goal is: Given an universe, Λ, of HMMs to chose the member, λ∗, which maximizes ∏ k=1..K P (Sk|λ), where S1, ..., Sk are the sequences we try to align. The procedure of choosing λ∗ can be thought of as training a HMM using training set {S1, ..., SK}. This is analogous to neural networks. The hope is that, after the training, the model have captured the essence of “what it means to be a member of this particular protein family”. If so, λ∗ will preferably generate amino acid sequences which belong to the family, i.e. sequences that exhibit the same